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Abstract 

To obtain the basis for combining various many-body techniques to QED in 
a consistent manner, we investigate the theory of quantum electrodynamical 
self-consistent fields. The reserch interest was born mainly of the electronic 
structure theory, thus we consider of atomic and moleculer systems as our 
main subjects. But the formalism is more fundamental. First, we derive the 
quantum electrodynamical Hartree-Fock theory by using the Green's function 
method. Then we construct a relativistic Hamiltonian written by creation- 
annihilation operators for electron and positron in a general form, and check 
that it reproduces the Hartree-Fock theory. We use this Hamiltonian to derive 
the time-dependent Hartree-Fock theory and random phase approximation, 
in the operator formalism. The relativistic Slater determinant of the Thouless 
parametrization is also used. Finally we discuss the applications and futher 
possible developments. 



*This work is firstly presented at the meeting of the physical society of Japan at Iwate univ., Sep. 
1999, and given as a part of the doctoral thesis of the author presented to Osaka univ., Nov. 2000. 
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I. INTRODUCTION 



It is natural to have an interest in the method to construct the relativistic many-body 
theory ( and relativistic self-consistent fields ). There are also many physical situations 
demanding the relativistic many-body theory for atoms, molecules and solids. Relativistic 
effects appear to be crucial in many cases. We have several relativistic extensions of many- 
body theory (T|-§[, but almost all of them are approximate theories from the beginning, 
not the fully relativistic formalism. The Hamiltonian is given under some assumption and 
approximation. Therefore they have validities and limitations. In the case of relativistic 
electronic structure theory, there are several approximation schemes. The famous Dirac-Fock 
method is based on the Dirac-Coulomb-Breit no-sea (the approximation completely neglects 
the effects of the Dirac sea) Hamiltonian |Ef. It has a Dirac type one-body Hamiltonian, 
and describes the interaction of electrons via the Coulomb potential and Breit operator. 
This Hamiltonian with several many-body schemes (Hartree-Fock, relativistic many-body 
perturbation theory, so on) give good numerical results for neutral heavy elements. But 
it can not be applied, for example, to describe the inner-core electrons of heavy elements 
or, electrons of highly ionized heavy elements such as lithium-like uranium. Because in 
these cases, quantum electrodynamical, QED, effects (such as the vacuum polarization or 
the electron self-energy) can not be negrected, we must take the Dirac sea into account. 
There is also a theoretical treatment for taking account of these QED effects ||. This 
method is based on the Furry picture of QED, and the QED corrections are estimated 
perturbatively with the renormalization procedure. This method gives very accurate results 
for highly charged ions. But it gives only correction. This method is not based on a single 
principle from beginning to the end: It must be combined with the Dirac-Coulomb-Breit 
no-sea scheme, for comparing with the experimantal value. Practically it treats only the 
one-particle QED effects, and the problem is reduced to that of hydrogen-like atoms. A 
matter of course, this method is only valid for the systems of a few electrons. 

A theory suitable for description of the interaction between electrons is intrinsically 
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QED |H . Electrons interact by exchaging a photon intermediately. Only QED can fully 
describe such a situation. Thus we have a fundamental question that how to construct a 
relativistic many-body theory based on QED in a consistent manner. And also how we 
can construct the self-consistent fields in QED, which it should be the central theme of 
the many-body theory in QED. In the nonrelativistic theory, there are various many-body 
techniques, and these methods are used in atoms, molecules, solids and nuclear systems. A 
large part of them is constructed in the operator method, and people use them in discussions, 
not only for stational properties but also for dynamical processes. It should be fruitful if 
these many-body techniques are applied to QED. It should enlarge the ability of QED, and 
should overcome the limitation of validity of it. By these considerations, we investigate a 
general formalism for the many-body theory in QED, especially the self-consistent mean- 
fields. In Sec.||, we derive the quantum electrodynamical Hartree-Fock (HF) theory by 
using the Green's function formalism. By this method, we can introduce the HF theory 
with no ambiguity. Thus we regard the result as the basis for our next step of investigation. 
In Sec.[]TJ, we give a relativistic Hamiltonian which is explicitly written by the creation- 
annihilation operators. We present a evidence that this Hamiltonian obtains the same result 
for the HF, with the Green's function method given in Sec.|T[ With this Hamiltonian, we 
examine the method introduced from the nonrelativistic many-body theory in the operator 
formalism, and derive the time-dependent Hartree-Fock (TDHF) theory , HF condition, and 
random phase approximation (RPA). We expect that the operator formalism will extract 
the capability of our theory. In this paper we mainly focus on the formalism. But we also 



discuss the possibilities of our theory, especially in applications of atomic physics, in Sec .[TV. 
In Sec.[V|, some general conclusion with summary, discussions of the remaining problems and 
the futher developments, are presented. We also present the finite-temperature quantum 
electrodynamical HF theory in appendix. 
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II. DERIVATION OF THE HARTREE-FOCK THEORY BY USING GREEN'S 

FUNCTION FORMALISM 



In this section we introduce a QED Lagrangian and Hamiltonian, and then we derive 
the Hartree-Fock theory by the field-theoretical Green's function method. Some parts of 
the contents were given in literature |5|-|9|]. But Eq.(21) is our original result. The form of 
Eq.(21) will become the basis of our theory developed in Sec.[TTI|. Here we intend to confirm 
the theory and to construct the basis for later investigation. We use the metric convention 
as ds 2 = g^ u dx^dx u = (dx ) 2 — (dx 1 ) 2 — (dx 2 ) 2 — (dx 3 ) 2 , and the metric tensors are defined 
as = diag(l, -1, -1, -1). 

The starting point for our discussion is the following QED Lagrangian 

C(x) = ~F^(x)F^(x) + ^(x)(tYd, - mo)^(x) + e^(x) 7 ^(x)(A M (x) + A^(x)), (1) 

where ip(x) designates the electron-positron field, m as its bare mass, and 7 M is the usual 
Dirac gamma matrix. The convention e = |e] is used. The anti-symmetric electromagnetic 
tensor is defined by F^ u = d^A v — d v A^. is the electromagnetic potential and A^ is the 
external field. To consider the atomic situation, we introduce the static Coulomb potential 
of nucleus and write it as Aq = U(r), where it is in fact nonrelativistic and classical. Now 
the coordinates describing the system are fixed at a nucleus which is located at the origin. 
As long as we treat the interaction between the nucleus and electrons by the external field 
approximation, it is difficult to describe the retardation effects. We suppose it is difficult to 
construct a fully Lorentz covariant theory in our situation. But it is enough for our aim in 
this paper. Here we also negrect the recoil effect of the nucleus. 

Adopting the action principle, we obtain the equations of motions (in the Lorentz gauge) 

d,F^ = (d 2 - V 2 )^ = -e^7"V, (2) 
{i-fdp - m - loU(r))^ = -erfA^. (3) 

Now we use the well known canonical quantization scheme for the Dirac fields 
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{ip a (x , x) , tpjj (x , y ) } = 5 a(S 5 (x - y ) , (4) 
{i> a (x , x), i>p(xo, y)} = {^i(x , x),$p(x , y)} = 0, (5) 

where a and (3 denote the spinor indices. We also quantize the electromagnetic field A^, but 
here we do not discuss the problem of quantization scheme with gauge fixing. For our aim, 
it is unnecessary. The equation of motion for photon are formally solved as 

A^(x) = -e d 4 yD^ v (x - y)${y)Yi>{y), (6) 



where D^x — y) is a full 2-point function for photon. Here we assume there are only virtual 
photons in the system. Next, we remove the photon degree of freedom and the equation of 
motion for the Dirac field becomes 

i 7 <9 ^(x) = {-il ■ V + m + lQ U{r))${x) + e 2 J d 4 yY^(x)D^(x - y)$(y)Y$(y)- (7) 

Then we find a Hamiltonian operator for describing the interacting many-electron-positron 
system 

H = J d 3 xtjj( y x)( y —i'j ■ V + m + j U(r))4>(x) 

+\e 2 J d 3 x J d^x^^D^x - y)i>(y)Y4>(y). (8) 

After these preparations, we derive the HF theory by using the Green's functions and 
Feynman rules. This method is famous in the textbook for the nonrelativistic many-body 
theory [|Kj, and here we extend it straightforwardly. In the perturbation theory based on 



the Dyson equations, the HF approximation is defined as an infinite-order summation of 
specific diagrams in expanding the fermion 2-point function ||-||[10]. Here we also obey this 



ansatz and take into account two diagrams as the lowest order vacuum polarization and 
fermion self-energy, in the central field. Thus we obtain the next formula for the HF total 
energy 

E = -iJ d 3 x(-ii • V + m + j U(r))tr(S% F (x, x)) 

-\e 2 J d 3 x J d'yDP(x-y),M^S^ F (x,x)MYS^ F (y,y)) 

+\e 2 J d 3 x J d'yDf{x - y^Ml^ F (x, y)YS^ F (y, x)), (9) 



where the first line is the kinetic energy while the second and the third are the direct and the 
exchage HF energy, respectively. The trace is taken over the spinor indices. Dp\x — y)^ is 
the free photon Feynman propagator, and especially in the Feynman gauge, it is written as 

iD F (x-yU = -ig^J (2?r)4 fc2 + v (10) 

Sp F (x,y) is the one-body fermion propagator constructed by the variationally optimized 
HF one-body functions. The fermion field operators are expanded generally as 



= E(^ +) (x)e-" x °a n + ^ } (x)e + ^ x ^ n ), (11) 

n 

k*) = E(& +) (x)e +ie " )xo ai + ^-)(x)e-^~H„), (12) 



where a n (ajj) is the electron annihilation(creation) operator while &„(&JJ is the positron 
annihilation(creation) operator. Here n indicates some quantum numbers suitable for the 
description of the system. Especially in the HF theory, ip^ becomes the HF one-body 
function and is its eigenvalue. The vacuum is determined by the relations 

a„|o) = o, (13) 

6 n |0>=0. (14) 
Hereby the HF fermion propagator is defined by 

tS F TF (x,y) = (%\T^ HF (x)^ (y)|$ > 



E^i +) (x)^ +) (y)e" le " (x °- yo) 0(x - y )e(e n - e F ) 

n 

n 

-EC ) (x)C ) (y)e +le " (xo - w ^(l/o-a:o), (15) 



where |<3> ) is the HF wavefunction occupied by electrons up to the Fermi level e F ( so called 
single Slater determinant ). Perform the Fourier transformation in time, it becomes 

sj-(x, y, e) = E^ +) (x)e } (y) f 0{e \ ~ €f) + 0{eF ~ en \ ) 

n \6 — 6n + IT) 6 — 6n — IT] J 

+ E 1>P (y) ■ ( 16 ) 

n 6 + 6n — IT] 



In (9), we also perform the Fourier transformation in time and employ the energy contour 
integration. For this purpose we must decide the integration path [[/j], so we will examine 
the structure of the 4-current because in (8) and (9), the interaction terms are written in 
the form of current-current interactions. If we choose the current operator which is normal- 
ordered about the Dirac vacuum, its vacuum expectation value vanishes in the zeroth-order. 
Consider, for example when we introduce an external field generated by a nucleus and 



describe the system by using the Furry picture 0,|TT|, the vacuum polarization is observed 
around the nucleus, so the vacuum expectation value of the current becomes non-zero even 
in the zeroth-order This situation becomes more clear if we concern to the fact that the 
zeroth-order in the Furry picture already involves an infinite order of summation in which 
each term is expanded by the external field. We want to include the vacuum polarization 
effect in our theory. For this purpose, the current generally have to be calculated by the 



Schwinger's expression [8,12 



2 \y-^x(x >y ) y-*x(x <y Q y (x-y) 2 >0 

The zeroth-order vacuum expectation value in the Furry picture becomes 



3a\x) = i-( lim + lim "]ti('-y^Sp\x, y)) 



2 \-+x(x >yo) y~*x{x <y o y 

\i<— m n>—m ' 

where iSp \x,y) = (0|T^(a;)^(y)|0). In this expression, in the extarnal field free case two 
sums cancel, but for the field present case two sums do not cancel completely but diverge (it 
must be mentioned that if there is no external field and the system is homogeneous like the 
electron gas, the normal-ordered current gives same result with above definition). Hereby 
this result should be interpreted as describing the vacuum polarization effect. To obtain a 
finite value, we must perform the renormalization procedure. Now the HF current with the 
vacuum polarization becomes 

if» = i?( 1™ + !/ m >r(7^# F (x,y)) 

l ^y^x(x >yo) y—>x(x <yo) / 



= -e E $ +) (*)7^i +) (*) 

-m<n<F 



e 
2 



E ^h^H*) - E & +) (*)7^i +) (*) ■ (19) 

\n<—m n>—m / 

With these observations, we conclude that the HF total energy should be written by the 
HF currents, given by (19). For this purpose we must treat the convergence factor, which 
reflects time-ordering of field operators in the propagator, to be symmetric and use both the 
upper and lower halves of the integration contour. Using the following equation 

E=-i J d 3 x J ^i(e ler h^+e- ier lim)(-i 7 - V + m + 7 t/(r))tr(^ F (x, z, e)) 
-ie 2 /^/^/$i/^4 0) (x-y,0U 



2vr J 2tt 



x I (e ieiT lim +e~ leiT lim)tr( 7 ^# F (x, z, e x )) 
xl(e ie2T lim +e- ie2T lim W-fSf F (y, z', e 2 )) 



+ ie 2 / rf 3 . / A / ^ / ^4 0) (x - y, Cl - e 2 )„ 

xtr( 7 ^f F (x, y, ^(e^ + e^ T )-fSf F (y, x, e 2 )^(e^' + e— ')), (20) 
and performing the energy integration, we obtain ( in the Feynman gauge ) 

E = ^J d 3 x hm (-«7 • V + m + loU{r) ) 

xtr( £ Vi" ) (x)^- ) (z)+ E ^ +) (x)e } (z)- E ^ +) (x)^i +) (z 

-e£~)<-Tri 

^ j <i\r j <f><) ^ 



-e£ } <-m -m<el +) <e f el +) >£F 



+ 



2 J J ^47r|x-y| 



x 

2 



\{^{ E ^- ) (x)^- ) (x)+ E ^i +) (x)^ +) (x)- E ^ +) (x)Vi +) (x))} 

xi{tr 7 ^( E W" ) (y)^i ( " ) (y)+ E W +) (y)$ +) (y) - E W%M (+) (y))| 



-m<e^'<e F e T '>e F 



1 2 /" ,3 /" ,3 expfileW -e, (±) ||x-y|) 1 f 



2 J J ™>™ 47r|x-y| 

X7^( E ^- } (x)^(y)+ E ^i +) (x)^ +) (y)- E ^i +) (x)^ +) (y)) 
xt"( E # } (y)^ } (x)+ E W +) (y)# } (x)- E W +) (yM (+) (x))) 



In (21), we handle not only iV-electrons but infinite number of interacting electrons and 
positrons. It is essentially different from the non- relativist ic theories. The direct term and 
the exchange term give the quantum correction self-consistently. In the direct term, there is 
no energy-flow and the photon propagator becomes an instantaneous Coulomb interaction, 
while for the exchange term the exponent in the propagator has the energy difference and it 
describes the retardation effect. Both the Dirac sea and the Fermi sea contain the Hartree- 
Fock quasi-particles. 

The (21) involves two types of divergences. One is the zero-point energy of the fermion 
system, while other has its origin in the interaction between particles. In principle, the 
zero-point energy can be removed by a simple subtraction, E HF = E — E Z p, and it gives 
the redefinition of the HF total energy. Here we use the next definition for the zero-point 
energy 

E ZP = -i [ d 3 x ( —-(e ieT lim +e~ ieT lim)(-z 7 • V + m )tr^ 0) (x, z, e), (22) 
J J 2tt 2 z ^ x z ^ x 

where iSp\x,y) = (0|T^(x)^(?/)|0) is the free, and no-external-field fermion propagator. 
On the other hand, the divergences which originate from the interaction must be treated by 
the renormalization methods in the quantum field theory [||]. For this purpose we should 
introduce some reasonable cut-off, and prepare the renormalization constants. But in the 
above case we can not apply the commonly used method for the field theory. There is a 
large number of efforts to perform the renormalization in atomic QED calculations 0. All 
of them are based on the perturbative QED. We have not met any of non-perturbative 
renormalization in QED, which can be applied to our thoery. This problem still remains as 
an open question. 
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III. OPERATOR FORMALISM 



The most part of QED in literature treats the translation and Lorentz invariant vac- 
uum, and which contains no external, central field. It is constructed by the field operators 
and Green's functions. On the other hand, the many-body theory is constructed not only 
by the field operators and Green's functions, but also by the creation-annihilation opera- 
tors, especially in the atomic and moleculer physics. In this formalism, we frequently use 
the occupation numbers for the one-particle states or the density matrices, and also the 
canonical transformation among the creation-annihilation operators take many active parts. 
These operators and the Hamiltonian are in the Schrodinger representation. We can find a 
large number of approximation schemes in this formalism, and they have various advanta- 
geous features. To combine these methods to QED, we investigate the creation- annihilation 
operator formalism for QED in this section. 



A. Hamiltonian 

This subsection is devoted to introduce a fully normal-ordered relativistic Hamiltonian 
for our theory. We essentially need a normal-ordered Hamiltonian, to use operator methods 
developed in the nonrelativistic many-body theory. We start from (8). This Hamiltonian is 
written in the Heisenberg representation, containing no approximation. The field operators 
are now substituted into the expanded form given in (11) and (12). For the kinetic energy 
term we can simply take the normal-ordering, and neglecting a c-number. Also the time is 
chosen as xo = to become the Schrodinger representation. About the interaction term, 
the photon propagator is substituted by the free propagator and we perform the Fourier 
transformation in time yo, and then we transform to the Schrodinger representation at the 
time Xq = 0. We use the normal-ordered currents and take a product of them to express the 
interaction term. Once more we employ normal-ordering for this product, and we obtain 
16 fully normal-ordered operators with 8 one-body operators. We simply negrect the latter. 
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Then we obtain the next fully normal-ordered relativistic Hamiltonian 

H = H K + Hi, (23) 
Hk = £ T$H^j + aJSt + ftiO,- - b]k), (24) 



i,j,k,l 

+a\b]a\ai + a|S]a^SJ' + a\b]b k ai + a\b\b]b k + a^a^ + aj^^a,- 

+biajb k ai + b\a 3 ^ k + b]a\biO,i + a\b]b\b.i + b]b k bi(ii + StSfSfcSi), (25) 

where we designate matrix elements as 

7j± = | d 3 x# ) (x)(-i 7 • V + m + 7ot/(r))^ (±) (x), (26) 



(27) 



Here, k, I denote quantum numbers for one-particle states, and Ae = \e k — e t \ = je, — €j\. 
Hereafter we use this Hamiltonian in our theory. It must be noted that this Hamiltonian 
(Eq.(23)) is in the Schrodinger representation. This Hamiltonian contains operators that 
they do not conserve the particle number, but conserve the total charge of the system, as 
the usual relativistic case. 

Next, we try to evaluate the expectation value of (23) for the Fermi sea, which is defined 

as 

l^)= II 4|0>- (28) 

— m<n<F 

It is clear that the naive application of usual nonrelativistic method, such as performing the 
factrization and taking the pairing, does not give the expression given in (21). Then here we 
propose to examine the relation between the contraction and the pairing. The contraction 
and pairing are defined by 

A'{t)B'(t') = T(A(t)B(t')) - N(A(t)B(t')), (29) 
A*(t)B*(t') = A(t)B(t') - N(A(t)B(t')), (30) 
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here • denotes a contraction, * denotes a pairing. The normal-order N is refered to the 
Fermi sea. The relation between the contraction and the pairings is 

A'(t)B'(t') = A*{t)B*{t')6{t - t') - B*(t')A*(t)9{t' - t). (31) 

The definition of the pairing gives (F\A(t)B(t')\F) = A*(t)B*(t'), the contraction becomes 

A'(t)B'(t') = (F\A(t)B(t')\F)6(t - t') - {F\B(t')A(t)\F}6(t' - t). (32) 

Especially in the case of t — t' 

(F\T(A(t)B(t'))\F) t=t , = (F\A-(t)B'(t') + N(A(t)B(t'))\F) t=t , 

= (F\A(t)B(t')\F)e(t - t')\ t=t , - (F\B(t')A(t)\F)9(t' - t)\ t=t , 

(33) 

Now we introduce the definition for the step function at the same times as 

6(t-t') t=t ,=9(t'-t) t=t , = ± (34) 

By this definition we obtain the following relation: 

(F\T(A(t)B(t'))\F) t=t , = (F\A-(t)B'(t')\F) t=t , 

= \(F\A{t)B(t)\F) - \{F\B(t)A(t)\F). (35) 

With the examination given above, we introduce a hypothesis that the Hamiltonian given 
in (23) is written by operators in the same-time T-products, and when we factorize the 
vacuum expectation value for an operator product with the aid of the Wick theorem, each 
contraction should be calculated by the definition given in (35). 
By using this assumption, for example 

Y t #% = Y t \{(M*i\F)-(F\*Aw) = \( E -E)*«> (36) 

ij ij —m<i<F i>F 

E-^I = Y.-\((F\b]k\F) - (F\kb]\F)) = \ £ (37) 

ij ij i<—m 
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ijkl ijkl 

= ~{ E -EW E -Efe (38) 

-m<i<F i>F -m<k<F k>F 

and other terms are also calculated. The final result becomes 

(F\Hk\F) = \ £ T~ + I( ^ -E)^t + , (39) 

i<—m -~m<i<F i>F 

(f\Hj\f) = Ui E -E)( E -E)(^t ++ -^ ++ ) 

L -rn<i<F i>F -m<k<F k>F 

+( E -E)( E )(v^~-^- + ) 

-m<i<F i>F k<-m 

+( E -EXEXW*-^-) 

-m<k<F k>F i<~m 

+( E )( E )(V£r" - vss~)}. (40) 

i<—m k<~m 

Now we can introduce the HF equations in the following form: 

+ E -E)(^E ±± -var + ) = +e{ +) , (41) 

k<~m —m<k<F k>F 

T u - + kY.+ E -E)(^ ±± -^ ± ') = -4' ) - (42) 

k<-m -m<k<F k>F 

These equations are rederived later, where we discuss the HF condition, and we can confirm 
that the definition of the HF equation given above has no conflict with it. The HF Hamil- 
tonian can easily be determined by use of the above results. By using (41) and (42), we 
obtain the HF total energy 

i<-m -m<i<F i>F 

4(E+ E -EKE+ E -EX^-^T") 

j<-m -m<i<F i>F j<-m -m<j<F j>F 

-E ZP . (43) 

Here we give the zero-point energy as in the (22), and it can also be calculated by the same 
way in this subsection. 
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Eq.(43) gives a consistent result with that derived by the Green's function method in 
the last section. The results given above imply that we can generalize this method to 
calculate expectation values for arbitrary products of operators. We always factorize the 
matrix elements of an operator product by using the Wick theorem into the sum of the 
products of contractions, and the contraction is calculated after the definition given above. 
It is clear from the above discussion that only the definition of the contraction for same-time 
operators is modified. Later we will use this definition to derive the HF, TDHF and RPA. 
It must be taken into account the fact that the operator method based on the Schrodinger 
representation, i.e., all operators appearing in the theory obey the same-time, t — 0. 

Now we assert that we obtain a basis to construct the relativistic many-body theory in 
the operator formalism. The Schrodinger equation in our theory becomes 

l*W) = EE^(*^--.^). ( 44 ) 

N 1 

i^Mt))=HMt)), (45) 

where is an iV-particle Slater determinant, and Af (t) is a time- dependent coefficient. 
It is clear from the above derivation, the Hamiltonian and the wavefunction depend on the 
gauge. A problem, how we should be treat the gauge dependence of our theory in many-body 
procedure, still remains. 



B. The Slater determinant 

In the Hartree-Fock theory, the many-particle wavefunction is expressed by a single Slater 
determinant, and the one-particle functions are optimized variationally. It is well-known that 
any Slater determinant can be written in the Thouless form 



|$) = e^|$ ), (46) 
G(z) = E^cfe, (47) 



where |$o) is a fixed A^-particle determinantal state. The parameter have to satisfies the 
Hermitian condition, = z* { , as e lS will be an unitary operator. The exponent generates 
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a canonical transformation on fermion operators. This second quantized representation of 
a single Slater determinant is very useful in operator methods. We try to generalize it to a 
relativistic form, because it plays essential roles in the later discussions. 
Let us introduce several Dirac vacua 

a n \0) = b n \0) = 0, (48) 
W = W = 0, (49) 
a^ F \0 HF ) = b* F \0 HF ) = 0, (50) 

where the first, second and third is old, new and HF vacuum, respectively. Determinants 
are given as 

l*o> = IR|0>, (51) 

n 

I^oHIRIO')- (52) 

n 

We want to obtain a relativistic exponential operator which transforms one to another 
representations. We assume the following relations for the canonical transformation: 





6 OjyyiC 


-iS 


n 


(53) 


4 


= e iS ale~ 


-iS 


n 


(54) 




= e*ble 


-iS 


n 


(55) 




= e iS b m e 


-iS 


= y^siPmnAn + 7mnM) 


(56) 



n 



where = S will be satisfied. Here new operators are expanded by a complete set of the 
old representation. The definition of the new vacuum becomes 

a' n \0'} = (e i§ a n e- iS )e i§ \0) = 0. (57) 

Therefore, the new vacuum is formally given by a superposition of the states in the old 
representaion. In the relativistic theory, the vacuum also be transformed. The relation of 
two determinants becomes 
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l*o) = 11^10') = 11(^4^)^10) = e^n^iO) = ^|$o). (58) 

n n n 

Now we introduce the exponential operator in the following form 

mn 

Here we demand that the matrices formed by the paramerters a mn , j3 mn , /3^ n , 7^ n should 
be Hermitian. Then the exponential operator formally given above is unitary It is clear 
that the (59) satisfies relations given by (53)~(56). Each parameter should be determined 
in the approximation scheme, and in the HF theory they are variationally determined and 
reflect the variational condition. Variations for one-particle functions are generated by each 
operator product in the exponent. It must be noticed that the second and the third terms in 
the exponent can be understood as the generalized Bogoliubov transformation to the Dirac 
vacuum, which does not conserve the particle numbers, and they describe the fluctuations 
of the Dirac vacuum. In the HF theory, the parameters included in them would be satisfied 
not only the HF variational condition but also vanish the matrix elements for the parts 
of the Hamiltonian which do not conserve particle numbers. These conditions must be 
satisfied for the diagonalization of the Hamiltonian. In (53)~(56), the new operators can be 
interpreted as quasi-particles, and they are given by superpositions of old operators. Later 
we use this exponential operator to derive the TDHF theory, the HF condition and the 
RPA equations. The most important point in this exponential operator is that it gives the 
degrees of freedom of a system explicitly. By this feature we can introduce many flexible 
approximation scheme, and can describe the dynamics of a system in detail, much more 
than the usual QED treatment. 

Now we must mention some mathematical points characteristic in the field theory. It is 
well-known fact that, in the theory of infinite degrees of freedom as quantum field theory, a 
different Fock space construct a separated Hilbert space, and two representations are uni- 



tarily inequivalent to each other flnj . Thus expansions used in (57) and (58) give just formal 
definitions. A vector of one representation can not be expanded by bases of another repre- 
sentation. For example, a state for the observed electrons can not be expanded meaningfully 
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by the bases of the bare electron Hilbert space. Thus we can not treat the interaction effect 
by only one Hilbert space. In fact, the renormalization procedure is a transformation to 
another representation. The renormalization method projects bare parameters to physical 
one, and transforms to the Fock space of observed particles. Therefore we interpret that a 
state of a new representation (the HF state) is formally expanded by old bases by using e tS , 
then we evaluate matrix elements, and have to employ the renormalization to transform to 
an observed particle (Hartree-Fock quasi-particle) space. In our theory, the Fock space for 
HF states is self-consistently selected with the renormalization procedure. 

Using the exponential operator given in (59), we can write the expectation value of the 
Hamiltonian (Eq.(23)) in a general Slater determinant (Eq.(58)) into an expanded form, by 
the same way as nonrelativistic cases. Here we will write a relativistic Slater determinat as 

|$(o:,/3,/r,7*)) = e^'^l^o)- (60) 
The expectation value of our Hamiltonian is given as follows: 

($(a,P,P*,Y)\H\<S>(a,P,P*,Y)) = (%\H\%) 

+i{$o\[H,S{a,P,(r,>y*Wo) 

+^(®o\[[H, S(a, (3, (3% 7 *)], S(a, (3, (3*, 7 *)]|<f>o) 

+0{S 3 ). (61) 

In (61), the first line is the HF total energy, the second line gives the first derivatives with 
respect to the parameters (a, /3, f3*, 7*), and it must be zero in the HF condition. The 
third line corresponds to the second derivatives and they determine the stability of the HF 
state, it is equivalent to the stability of collective modes in the RPA. It is obvious fact that 
the operator formalism makes possible to do these discussions. By using our theory, for 
example, it may be possible to discuss the instability of the HF solution, or to apply the 
coupled-cluster method in the electronic structure theory. 
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C. The time-dependent Hartree-Fock theory 

In this and in the following subsection, we will derive the TDHF, HF condition and RPA, 
by using the operator method, based on the previous discussions. Here it is also the aim 
to demonstrate the usefulness of our theory. We take the starting point to the TDHF, the 
theory for the time evolution of the mean-field [ |T5|j . 

Using the Hamiltonian ( Eq.(23)), we construct the action functional 

= - iJ|«>(c(t), 7 *(t))> (62) 

Here the time-dependence of the determinantal state is borne on the parameter set 
(a, /3,/3*, 7*). Introduce the Dirac-Frenkel variation principle, 

SShf = sj*' dt(®{a{t),P{t),P*{t) n *{t))\^ ~ Hma(t),/3(t),P*(t),y*(t))) = 0. 

(63) 

Thus we obtain the stationary condition 

Bi 9 

= 5ma,(3,^,j*)\i— - tf|$(a, /?,/r, 7 *)>. (64) 

Therefore we take the derivatives with respect to the parameters and derive a set of TDHF 
equations 

Bi 9 

= ($(«, /3, /T , 7 *) | [ala n , i— - H]\$(a, P, /3*, 7*)), (65) 

Bi 9 

= ($(«, /3, P\ 7*) I [SU>1 i-Qf ~ H }\^ P, /?*> 7*)), (66) 

Bi 9 

= ($(«, p, p*, 7*) I Kan, - mm<*, P, P*, 7*)), (67) 
~ ~ Bi 9 

= ($(«, p, P\ 7*) I Ml, i-Qf - #P(«, P, F, 7*)), (68) 



where 



BiS 1 • - *± 



mi; 



Another form can be obtained for the TDHF equations. We also derive the TDHF equations 
in the density matrix form because later we use them for the derivation of the HF condition 
and RPA equations. Hereafter we write |3>(/c)) for \&(a, (3, (3*, 7*)). The definition of the 



Pnm( K ) Pnm( K ) 
y Pnm( K ) Pnm( K ) j 



(70) 



density matrix is 

P (k) 
where 

i 2 * - 

= ($ 1 a% a' n I %) + i ($ 1 > S \ I $o> + g - <$o I W m a n , S] , S]\<S> + --- 

= P^ (0) +P^ (1) (-)+P^ (2) («) + •■-, (71) 

where an expansion at a specific point is given. The definition for others are, respectively, 

p^(«) = ($(«)|at l St|$(«)>, (72) 

/»nm(«) = <^(«)|6mOn|^(«)>, (73) 

p~(«) = ($(«)|Uil$(«)>. (74) 

It must be noted that a* a m , alp^, b n a m , bj)^ construct a Lie algebra, so it is not necessary 
to concern ourselves about other types of operator-pairs throughout our theory. Using the 



relations 



g~ t PuM = (m\[ala n , ^]|^(«)>, (75) 
($(«)|ff|$(«)> = ($(«)|[ff,at l a B ]|$(«)>, (76) 



d%OL mn 

and so on, we obtain the TDHF equations in the following forms: 



= (78) 

= -g^r(^)\HMK)), (79) 

'' m PnM- 
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J nm 



i^PnM = " or— ($(«)|ff|$(«)>. (80) 



In our formalism, theories of beyond-TDHF also can easily be constructed. For example, 
we can introduce a multi-configurational wavefunction \^(t)) = Y^i Aj(t)\§j) (here |$/) is a 
Slater determinant) for making an action functional and adopt the Dirac-Frenkel variation. 
Then we obtain the multi-configurational TDHF in QED 

5S MC = 5 J dt^AfflMtWA^^'^m ~ ^e^'^'^^j), (81) 
where the configuration coefficients and one-particle functions are optimized simultaneously. 



D. The HF condition and the random phase approximation 

The TDHF gives a picture that the presentation point for the system moves in a param- 
eter space. The RPA has the meaning of the linearized approximation to the TDHF. It is 
well-known that the RPA equation can be derived by several different methods. The famous 
equation-of-motion method developed by Sawada [TJJ can also be applied to our case, but 
in our case it is too cumbersome. Here we use the TDHF equations and derive the RPA 
equations by expanding the density matrices ]T3| . This method is convenient to derive the 



RPA and also we can easily interpret the physics and relations to other theories. 

First, we try to write down the expectation value of the Hamiltonian (Eq.(23)) in a single 
determinant |$o)- It is given as follows 

= ]T7;f { P ++(k) + />/f (*) + P ;+(k) + P ;r( K )} 

ij 

+ 2 X^(^j = W ±± _ Vitkj ±± ) 
ijkl 

x{pJ + («0 + + Pl + ^) + /<£•"(«)} 

xK + («0 + PiST(«) + PTk + (K) + Pi h 7(K)}. (82) 

Next, we introduce an exponential operator to use the expansion of the density matrix by 
parameters. This operator must explicitly divide the one-particle space into above and below 
the fermi level. Then we introduce the next form 
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e lS = expj if ^^VajA, + E a pg^g + E^w^jA + tfiptyflv) 

I ^ \u> pq HP 



+ E(/3pmaJ&m + Ppmbmap) + E^^m + ^nfi^n) + E ImrfiA ) \ ■ (83) 



""p"m 1 ypm^m^pj i ^^(im l «( 1 l ' m 1 y H m m L 1 ' 
pm fxm ran 

Hereafter, we use /i, z/ for electrons above the Fermi level, p, q for electrons below the Fermi 
level, and m, n for positrons, and i, j, k, I denote arbitrary states. In the above equation, only 
the third and the fifth terms in the exponent give physical effects, and the others describe 
only the phase degree of freedom. Then we negrect the latter and 

expli(j2(a up ala p + a* p aja^) + ^(^ m ft + P* m b m aM L (84) 

(, ^ HP um ' ) 

We will use this operator. Now we expand the density matrix in (82) by the parameters, 
and rearrange them in the order of the parameters. As the Oth order terms we obtain 



ij ijkl 



and it corresponds to ($ \H\$ ), the HF total energy, given in (39) and (40). Next the 1st 
order terms are given as 



E^f{pJ +(1) (*0 + pt {1) M + p7, +(1) (*0 + p7 w ^)} 

ij 

ijkl 

x{pj +(1) (-) + pr (1) (-) + p. +(1) (*o + pv w ^{pf + P7 k - { % 



(86) 



which correspond to i(& \[H, S]\& ). Calculate each matrix element and rearrange them, 
and we obtain 

y~X-^V ^ a HP ~ iQ> up ) + {T m piPiim — T a m ^finm) 

HP n n 

U E+ E -E)(e((^-^)<« w -(^S^ 



2 



K-m -m<l<F 1>F / \ HP 



-VX ±+ )Kp) + E^nS" - V ^UH ±+ )^m - (V;-^ - V+±±-)i(3; m ) . (87) 



Hm 
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Here if we take the derivative with respect to parameters, then we obtain the HF conditions 
(they are the same as the HF equations) and which coincide with previous definitions given 
by (41) and (42). The second order terms become 



£7Tor (2) («) + p;r (2) (*o + ^ +(2) («) + p;r (2) (-)} 

ij 

ijkl 

x (pJ +(2) («) + («) + («) + %" (2) («)) 

x(riT (1) («) +rir (1) («) +ft fc _(1) («))}■ (ss) 



■ 2 A, A, /\ 

It corresponds to |y($ | [[-ff, S], S]\<& ). We also calculate each terms, take into account 
the symmetries of matrix elements, Hermiticity of parameters, and choose the diagonal 
representation for the HF Hamiltonian, then we obtain 



E (^wt + + o ( E + E " £)0 / wii« ±± ^/i«t ±+ ) 



9 

MP i<-m -m<l<F 1>F 

-t p + p + - \{ E + E - E)(^T ±± - ^ ±±+ ))«; P ^ 

Z<-m -m<l<F 1>F 



+E( T - + + 9 1 ( E + E - EXC^ - K^ ±+ ) 

» Z<-m -m<l<F 1>F 

~T m m ~ E + E ~~ E)(Knmzf ± — KnZlm ))fltmflvrn 
K-m -m<l<F 1>F 

E E{(^/W ~ Vfivqp ) a np a vq 

~c^Yp[xqv Vpvqfi )(^fJ,p(^iyq l^Y [ipvq ^fxqvp ) a fip a uq } 

- /3* — —TV" — 1 — h — V" — ^ — h ")/3 8 \ 

2 \ /imra '[mum irptmyvn ^\"mfmu ravnp, jy^rnyunj 

+ W(fV +++_ -V r+_++ )a 3* - (V +++ ~ - V + ~ ++ )a* 3* 

fip vn 

"K^rov + ~~ ^tump + ) a p,pfti>n ~ (Vp + ^v + ~ ^pun/x + ) a /JtpPvn\- (89) 

Now the TDHF equations from which the RPA equations are derived become 
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= -^<$(«)|^|$(«)>, (90) 

MP 

4^(«) = -»J-^(«)l^l^(«)>. ( 92 ) 



= -o7^W«)l^l*(«)>- ( 93 ) 



dt rmvv ' dip, 

■ d_ + _, s 9 

'•dt PvmKK) ~ dtp 

Expand both sides with respect to the parameters as 

4KV" (0) +^ +(1) w+p;; (2) («)+•••} 



j($o|# |$ > + i<$ |[jf , S]|$ > + l y (®o\[[H, S], S]\%) + •••}, (94) 



and so forth. In this expansion, the first term of each side is a constant, and the second 
term in the right-hand side becomes zero when HF condition is satisfied. Therefore, taking 
the lowest order non-zero matrix elements of both sides, we obtain 

<,,(k) = -«y «* < 95 > 

We assume that the system is in small normal-oscillations at a stationary point, and the 
number of the normal-oscillations equals the number of one-particle excitations. Then we 
can expand each freedom by the normal-oscillations 

« w = E(^ e ^ fct + Gi e+ *' fct )» (99) 

k 

< P = E(0 + ^ + Q>~ lWkt )i (ioo) 
k 

(3 vm = J2(Rt m e-^ + S k um e +i ^), (101) 

k 

film = E( i & e+iWfct + St* m e-^ kt ). (102) 
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Adopting them to the former relations we have derived, and picing up by the same frequency, 
we obtain the following relations (with complex conjugates) from (95)~(98) 



j->k 


— f e (+) _ e W)p ft 
\ /J, p > up 














. V^f (\r++++ 
' 2-^i \ v fipqv 

vq 


t /++++\ r>k 
fivqp ) vq 


/T/ + + + + 

V * fipvq 


T ^ + + + 4 
fiqvp 


)Qvq) 






^ ; { ( ^fipnv 
vn 


T/+H h\ T>k 

v uvnp ) ^vn 


- (v +++ - 

\ ' upvn 


- v + ~ + - 

' unvp 


f \ nk* 1 
!°vn\i 


(103) 


-UkQ% 


~ \^u ~ e p )^ up 














' Z_^i \ \ V fipqv 
vq 


uvqp l^ivq 


— (V ++++ 

\ " fipvq 


— 1/ 

fiqvp 


Y\ jyk* \ 
) P uq) 






^ ^\^fipnv 
vn 


y++-+)S k 

uvnp J vn 


- (V +++ ~ 

V ' upvn 


unvp 


) K vnjy 
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ok* 


I _(+) i A—)\ ok* 
\ v ' m )um 














+y\(v- ++ - - 

' / j [\ Y mvun 
fin 


y—++)S k * 

Y mnfiv i^fin 


- (v- + - + 

\ ' mvnfi 


- v- + -^ 

Y mfinv 


~)R k ] 

1 1 fin J 






+ Ti(v~ +++ - 

1 / j [\' mvfip 
up 


V~ +++ )O h * 

mpuv / up 


- (V~ +++ 

\ mvpu 


- v~ ++ - 

mfipv 


h\pk 1 
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= f e (+) _L e (-))p fc * 
\^v 1 / Vm 














+ Y\(v~ ++ - - 

' / j \ \ y mvfin 
fin 


V R k * 

" mnfiv ) ^fin 


- (V~ + ~+ 

\ mvnfi 


- - V~ + ~' 

' mfinv 


f )S 1 fin\ 






+y\(v ++ - + - 

1 / j y V ' upmv 
up 


T/+H V\ pk* 

fivmp ) flip 


- (v ++ - + 

\ pfimv 


- v ++ - 




(106) 



These are the general eigenvalue equations in our RPA theory. It is easy to obtain the 
polarization function 

nJSfry.u;) = -ij ^tr[ 7 ^f F (x,y,e + o;)7^f F (y,x,6)], (107) 

from our RPA equations by use of a procedure employed in the nonrelativistic theory. These 
equations treat the Dirac vacuum on an equal footing with the iV-electrons, as a many-body 
systems. 

Ionescu et al. also obtained an RPA equation for Dirac vacuum []17|| , by using the method 
of Bethe-Salpeter equation |I(J. They studied the collectivity of Dirac vacuum around heavy 
nucleus, but could not find the evidence of it. 
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IV. SUGGESTED APPLICATIONS 



In this section, we propose some applications of our theory. Because it is besed on a 
simple and natural formulation, we should find many physical situations for our theory. Here 
we limit our discussion to the atomic physical situations. 

In the context of the electronic structure of atoms, there are three effects in atoms: The 
electron correlation, the relativistic effect and the QED effect . The electronic structure of 
atoms is determined by the relation of these three factors. The electron correlation depends 
on the electron numbers and can be treated by several many-body methods. The relativistic 
effect becomes large with increasing the atomic number. This can be treated satisfactorily 
by the Dirac-Coulomb-Breit no-sea scheme. The case of highly ionized heavy atoms, the 
QED effect becomes dominant and it can not be neglected. This effect is treated by the 
perturbation theory in the Furry picture of QED. But practically this method can treat only 
the system of atoms with a few electrons. In fact, especially in heavy atoms, as the ionicity 
becomes high and the electron number decrease, the many-body effect becomes small and 
the QED effect becomes large. Therefore we argue that our theory should be applied to the 
cases where the electron number is not few and the many-body effects must be take into 
account at least in the mean-field level ( can not apply Furry picture perturbation scheme, in 
practice ) , also the QED effects are strong and can not be neglected. Thus heavy elements 
in middle level of ionicity should be one of the subject of our theory. The Dirac-Coulom- 
Breit no-sea scheme has been applied to only neutral or almost neutral atoms. The method 
of the QED correction works only for the highly ionized heavy atoms. Our theory should be 
suitable for intermediate case between them. We can make mention that only parts of the 
whole domain of the atomic physics was studied. 

To explain the excited states or optical phenomena of the system, we need the information 
about the many physical properties, such as the excitation energies, transition amplitudes, 



multipole polarizabilities and oscillator strengths ||19|| . The RPA can treat these properties. 
The RPA method contains specific diagrams in all orders. Thus generally the RPA gives high 
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numerical accuracy. The inner-core excitation or ionization in heavy elements like uranium, 
the QED effect is strong and we need a quantum electrodynamical theory. Our RPA ( Eqs. 
(103) ~ (106) ) should be applied. Our RPA is suitable for the phenomena where both the 
electron correlation and QED effect is important. 

In the heavy-ion collision, if the nuclear charges are large enough, a spontaneous electron- 
positron pair creation occurs. Our RPA can give all types of excitations ( electron excitation, 
ionization, pair-creation ) which will be happened in the heavy-ion collision. Futhermore, in 
these collisions, the adiabaticity breaks down. In such a strange chemical reaction we should 
treat the dynamics of the electrons and nuclei on an equal footing. There is a method that 
does not use the adiabatic potential energy sufaces but uses the time-dependent variational 
principle [[HJ. Our TDHF ( Eqs. (77)~(80) ) can be combined with such a method to 
describe the heavy-ion collisions. 

The problem of the renormalization still remains. But if it is solved in the future, 
our method should gain a large ability for numerical calculations, and should find various 
applications especially in the atomic physics. We should obtain several maps or charts for 
various physical quantities in the atomic physics, which they have two coordinates ( atomic 
number and ionicity ). 

V. CONCLUDING REMARKS 

In this paper, we have investigated a method to construct a quantum electrodynamical 
many-body theory, especially the mean-field theory. We have given the form of the HF in 
QED, and also some problems were discussed. We have taken the result obtained there as 
a starting point of our study, and have investigated the operator method. A relativistic 
Hamiltonian was derived on QED. The relativistic Slater determinant of Thouless form 
was obtained. We have discussed that in our case, the Wick theorem should be modified. 
After those preparations, we have introduced the TDHF and RPA in our theory. In our 
observation, the theory showed no difficulty at least in the mean-field level. All of our theories 
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include the Dirac sea as the dynamical degrees of freedom. This point is the characteristic 
feature of our formalism. This feature extends the capabilities of our theories than others. 

We should point out the possibility of the path-integral formalism. In a non-relativistic 
theory, Kuratsuji and Suzuki developed the method using generalized coherent state (GCS) 
representation JHJ]. The crutial point of this method is using GCS for obtaining the resolu- 
tion of unity and construct a kernel in path-integral form. They gave a proof that the Slater 
determinant of Thouless form can be used as GCS. We propose to apply this method to our 
theory. In the imaginary time formalism, the partition function becomes 

Z = Tre~ pA 

= J D(a,/3,/3*,7*(r))exp[- jf" dr(a, 0, /T, l*(r)\-^ + H\a, 0, /T, T *(r)>] , (108) 

here \a, (3, /3*, 7*) is the Slater determinant and r is the imaginary time. The path is de- 
termined by the movement of the coordinate (a, f3, (3*, 7*). This kernel is exact for our 
Hamiltonian. If we apply the stationaly phase approximation, we obtain the time- dependent 
veriation principle leading to the TDHF equation, (63). Expand in the vicinity of the HF 
state by parameters, the following form is obtained: 

Z = e-? E «F yp(a,/3,/3*, 7 *(r))exp[- £ dr{($ \[S, 3]|* ) - ~($ \[[H , S], S]\$ )} . 

(109) 

About futher progress, we mention about the method of beyond-RPA. The RPA describes 
the small oscillations at a stationary point and only suitable for describing the low excita- 
tions. We can investigate, in principle, the theory for high excitation and large amplitude 
anharmonic oscillations, such as the mode coupling theory. This theory should be applied to 
the system which is under both the strong QED effect and strong quantum fluctuation. We 
can also extend our theory to the non-equilibrium statistical mechanics through the TDHF 
formalism, the same way as the nonrelativistic theory. The works for these problems should 
become in the future. 
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APPENDIX: THE FINITE-TEMPERATURE HARTREE-FOCK THEORY 

In this appendix, we derive the quantum electrodynamical finite-temperature Hartree- 
Fock theory by using the temperature Green's function formalism. To our knowledge, no 
finite-temperature HF theory has been presented in literature. But of course, there is the 
finite-temperature QED in literature p2| . 

First, we give several definitions. The grand partition function and statistical operator 



are given as 

Z G = e -? n = Tre~ pk , (Al) 

p G = Z G l e-? k = e^~ k \ (A2) 

The Hamiltonian for the system is given as 

k = k + kj, (A3) 

K = H -fiN = J d 3 a^(xr)(-z 7 ■ V + m + lo U{r) - 7o /z)</>(xt), (A4) 

K I = H I = -eJ rf 3 x|(xr)7^i M (xr)^(xr). (A5) 

The field operators in the Heisenberg representation 

^(xr) = e* T 4(x)e"* T , (A6) 

V^(xr) = e^4(x)e-^, (A7) 

i^(xr) = e^i,(x)e-^, (A8) 



where the first, the second and the third is the fermion, its Dirac conjugate and photon, 
respectively. The fermion field operators in the interaction representation are defined as 
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= EW +) ( x )e" (£ ^" )T fl, + ^" ) (x)e + ^" ) +") T St), (A9) 



^(xr) = e^ OT ^(x)e- 75 " OT 



= E(^ (+) (x)e +{ ^ M)T a] + ty-\x)e-$ ^Hj). (A10) 
j 

Then the free fermion temperature Green's function becomes 

<S(°)(xt,x't') q/3 = -(T T [^(xr)^(xV)]) 

= -{«(r - /)E(f(x)f(x')e-^-»( T - T '»(l - 

3 

+^" ) (x)^- ) (x')e + ^" )+ ^( T - r 'V(4" ) )) 

(+)^,7;(+)^'U-(4 +) -M)(r-r')f^(+)^ 
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4-\x)$-\x')e + ^ )+ ^- T '\l - /(e^)))}. (All) 
Here each fermion distribution function becomes 

= , /(4 _) ) = rr~ , (A12) 

and (O)o denotes the ensemble average e /3Q °Tr{e~ /3i ^°0}. Because of the fermion statistics, 
S satisfies the anti-periodic condition. We perform Fourier transformation in imaginary 
time. The Matsubara frequency will be odd, u n = (2n + l)n/(3 and we obtain 

5 (0) (x,x',cu n ) = f We^S^fx^r) 
Jo 

= -E ^ +) (x)^ (+) (x) jf' dre^e'^-^il - f{e^)) 
- ]T V^WV^V) \l dre <WnT e +(e i" )+ ' t)T /(e5" ) ) 



_ E / ^ (+) (x)^(xQ | (Ai3) 

The free photon temperature Green's function in the Feynman gauge is given by 

T?(°)(xr,xVV = -(T r 4(xr)i i ,(xV)) 
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1 f d k ^k-(x-x')-ia; n (T--r / ) 



/^J (27T)3 W 2+k* 
1 -|w„||x-x'| 

9-4^1 e "'"" ,T "°- < A14 » 
where V satisfieds the periodic condition because of the Bose statistics, and frequency be- 
comes even u n = (2n)ir/(3. 

Now we give the HF internal energy on a basis of the finite-temperature Feynman rules 
and also of the results of previous sections. Take the symmetric form for the convergence 
factor, and we obtain 

Ehf{T, V, fi) = J d 3 x^^(e-^hm+e-^lim) 

x (-ry • V + m + -y U (r) - -fofi)trS HF (x, z, u n ) 
+ le 2 / d 3 x J d»yl £ I £ 2><°> (x - y, 0),, 

x 1 (e*""" lim +e-^" lim) tr( 7 ^ HF (x, z, w n )) 
x^fe*""'" lim + e - iw »"' lim)tr(7^ F (y,z',cu n 0) 

2 ^ z'^y z'— >y' V / 

xtr{7^ l S HF (x,y,^„)i(e^ + e~ iu> ^) 
xYS HF (y, x, c^O^e*"*'" + e"^'")} 

-E ZP , (A15) 

where <S^ F is the self-consistently determined HF temperature Green's function. Here the 
definition of zero-point energy is the same as (22). The internal energy depends temperature 
T and chemical potential //. Perform summations for each frequency by the usual complex 
energy plane method. We obtain 
1 r 



E HF (T, V,n) = -J d 3 x Jim (-17 • V + m + 7o C/(r) - lol i) 

xtr(E^ (+) (x)^ (+) (z)(2/( e 5 +) ) - 1) +E^" ) (x)^" ) (z)(l -2/(6 



(-) 
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x^r{ 7 ^(5:#)(x)#)(x)(2/( e ( + )) - 1) +Et H Wf'W(l -2/(6^)))} 
x^r{ 7 ^(^^ +) (y)^ (+) (y)(2/( e ( +) ) - 1) + E^y^yXl - 2/Ce^)))} 



x7^(E^ (+) ( x )^ (+) (y)(2/( e ! +) ) - 1) + E^ _) W^ H (y)(i - 2/( e 



exp(i 


6 (±) - 




|x-y|) 


47T 


x- y 





(-) 



x7^(E^ (+) (y)^ (+) (x)(2/( e 5 +) ) - 1) + £^ -) (y)$ -) (x)(i - 2/(4-^ 

-E ZP . (A16) 

The above equation includes the effect of the Dirac vacuum correctly. Take the zero- 
temperature limit, it recover the zero-temperature theory given in the previous sections. 
The Dirac vacuum also has thermal energy, depending on the temperature and be treated 
in the same footing as many-particles. (A16) also has the renormalization problem. Our 
result given in this appendix should be a starting point to the next step, for example, the 
operator formalism for finite-temperature mean-field dynamics. 
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